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Abstract 

A new result enables direct calculation of thermoelastic damping in vibrating elastic 
solids. The mechanism for energy loss is thermal diffusion caused by inhomogeneous 
deformation, flexure in thin plates. The general result is combined with the Kirch- 
hoff assumption to obtain a new equation for the flexural vibration of thin plates 
incorporating thermoelastic loss as a damping term. The thermal relaxation loss is 
inhomogeneous and depends upon the local state of vibrating flexure, specifically, the 
principal curvatures at a given point on the plate. Thermal loss is zero at points 
where the principal curvatures are equal and opposite, that is, saddle shaped or pure 
anticlastic deformation. Conversely, loss is maximum at points where the curvatures 
are equal, that is, synclastic or spherical flexure. The influence of modal curvature 
on the thermoelastic damping is described through a modal participation factor. The 
effect of transverse thermal diffusion on plane wave propagation is also examined. It 
is shown that transverse diffusion effects are always small provided the plate thickness 
is far greater than the thermal phonon mean free path, a requirement for the validity 
of the classical theory of heat transport. These results generalize Zener's theory of 
thermoelastic loss in beams and are useful in predicting mode widths in MEMS and 
NEMS oscillators. 

1 Introduction 

High Q resonators are central to the development of new devices and applications that 
include RF filters, next generation MRI systems, and torque magnetometers. Silicon based 
micro- and nano- electromechanical systems (MEMS/NEMS) oscillators are the candidates 
of choice, and include free standing planar devices, such as double paddle oscillators (DPOs), 
and micro-cantilevers. As the oscillators shrink in size, it has been found that the Q achieved 
is orders of magnitude smaller than expected based on classical fundamental loss mecha- 
nisms. Many mechanisms have been proposed, including surface loss m2IS| that increases 
with the surface to volume ratio. However, under controlled conditions with minimal surface 
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defects and adsorbates, measurements on silicon DPOs have shown that room temperature 
losses are adequately described by thermoelastic relaxation, while unexplained mechanisms 
operate at lower temperatures. Interestingly, the mode of vibration of DPOs is designed to 
be primarily torsional with very little flexure (and hence no thermoelastic coupling) . How- 
ever, as demonstrated by Photiadis and his co-workers ^ |31 El E] it is precisely the small 
amount of flexural motion that accounts for loss in these supposedly torsional oscillators. 

The purpose of this paper is to provide a consistent theory for predicting intrinsic dissi- 
pation arising from thermoelasticity in elastic structures. Particular attention will be given 
to flexural motion of thin plates. This work is a step towards understanding the fundamental 
limits of dissipation in small structures such as NEM and MEM devices. 

We are concerned with determining the thermomechanical loss of elastic modes, for 
example, the flexural mode of a rectangular thin plate. A useful point of departure is the 
classic theory of Zener 8 for anelastic thermoelastic damping. The key to this approach 
is the assumption, implicit in Zcner's work, that there is little relative difference between 
the isentropic (unrelaxed) and isothermal (relaxed) mechanical responses, and hence the 
mechanical and thermal problems are essentially decoupled. Since the thermoelasticity is 
weak, the transition from the instantaneous or unrelaxed system to the relaxed state can 
be viewed as a quasistatic thermal process, governed by the standard equations for thermal 
diffusion, although now in the presence of an inhomogeneous deformation. 

Energy loss in a mechanical oscillator is measured in terms of the quality factor, defined 
as Q = 27t£o/ A£, where £o is the mechanical energy of the oscillator and A£ is the energy 
loss per cycle. The quality factor for a lightly damped single degree of freedom system with 
nondimensional damping C <C 1 is Q = l/(2C)i ^tnd by assumption we only consider systems 
with light damping, or Q ^ 1. The relation between Q and aat, the attenuation per unit 
length of a propagating wave of frequency ui, is Q — uj / {2aatv), where v is the speed of 
energy propagation, also equivalent to the group velocity. This identity can be derived by 
assuming the energy is a quadratic function of the field variables, so that energy decays 
with distance d as e"^""*'*. The distance travelled in one period is d — 2tiv/ijj, and hence 
the fractional decrease in energy per period is 1 — e"*'^""'^/'^ « ATraatv/oj from which the 
relation for Q follows. 

Thermoelastic loss can be most simply viewed as a relaxation mechanism with a single 
relaxation time t. The generic frequency dependent quality factor Q{lo) for a relaxation 
mechanism is 

Co 1 -|- w^r^ 

where Ac is the (relatively) instantaneous increase in elastic modulus, cq cq -I- Ac, caused 
by the process under consideration. The change in elasticity is well known for thermoe- 
lasticity, and r has been estimated for several configurations. Thus, Zener [S] showed for 
flexure of a beam that 

Ea'^T ujTo 

Cp 1 + Uj'^T^ 

where E is the isothermal Young's modulus, T the absolute temperature, a the volume 
coefficient of thermal expansion, and Cp is the heat capacity at constant stress. The char- 
acteristic relaxation time is tq = h'^Cp/{TT^K) where h is the thickness and K the thermal 
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conductivity. In fact, as Zener first demonstrated [J, the simple expression (0) is the lead- 
ing term in an infinite series which is well approximated by the single term (see and 
Appendix B). Zener's method was derived in the context of scalar problems, where the 
strain, for instance, involves a single component. An important example is the flexure of 
a beam or reed, as considered originally by Zener, and later by others, for example, |l(Jj . 
The present work generalizes Zener's method to consider general elastic deformation. This 
includes the inhomogeneous deformation associated with modes in thin plates and other 
structures. 

Since the original work of Zener numerous papers have appeared on thermal relaxation 
in the context of coupled thermoelasticity. In a pair of papers m Alblas provided 
a rigorous formulation using continuum thermomechanics and linear elasticity theory for 
isotropic materials. He derived detailed and explicit solutions for the thermoelastic damp- 
ing in vibrating beams, including the circular rod and the rectangular beam. The result 
for the latter was derived separately by Lifshitz and Roukes although Alblas' solution 
is the more general of the two. These analyses are compared with the present formula- 
tion later (see Appendix B). Kinra and Milligan ^1] again derived the coupled isotropic 
thermoelastic equations and provided a solution for unidimensional structures, including a 
discontinuous beam. Perhaps the most thorough analysis of thermal damping in the con- 
text of the coupled equations of thermoelasticity is due to Chadwick By considering 
a modal decomposition of the elastic and thermal fields, an exact relation for the complex 
valued frequency of oscillation of each mode was obtained. This enabled Chadwick to derive 
a generalization of Zener's expression for the thermoelastic damping of an arbitrary elastic 
body. Chadwick subsequently derived the governing equations of thermoelasticity for thin 
plates and beams ^B]- The equations are in the form of coupled equations, one of which 
reduces to the classical equations for the structural mode, for example, flexural waves in 
thin plates, and the other the temperature diffusion equation, in the limit of zero coupling. 
The analysis in |15[ll6j is restricted to isotropic solids. 

This paper has several objectives. The first is to demonstrate how the Zener model 
follows from the full equations for the coupled dynamic system by using a consistent ap- 
proximation scheme. In the process we generalize Zener's approach to incorporate general 
elastic deformation, specifically the elastic stress and strain tensors. The main applications 
are to thin plate structures, for which we obtain a Zener-like result for arbitrary flexural 
deformation that includes the general curvature tensor. Our results will also include the 
possibility of thermal diffusion in the lateral direction in thin plates, which is explicitly 
ignored in Zener's approach, but was considered by Alblas ^T) and Chadwick 
However, it will be shown that circumstances under which lateral thermal flux becomes 
important coincide with the limit in which the thin plate theory is no longer applicable. 

The paper is arranged as follows. Governing equations of thermoelasticity are presented 
in Section 2 for anisotropic elastic bodies. General solutions are discussed in Section 3 with 
no particular type of structure in mind. The theory is applied to thin plates in flexure 
in Section 0] and a non-dimensional modal participation factor (MPF) is introduced which 
defines the local contribution to thermoelastic (TE) loss in terms of the plate curvature. 
An alternative method for deriving the TE loss of travelling flexural waves is presented 
in Section [S] using generalized plate equations. The effects of lateral thermal diffusion are 
discussed in the context of travelling wave solutions in Sectional 
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Table 1: Thcrmoelastic constitutive equations. The left column indicates the dependent 
variable. Each subsequent column corresponds to a set of constitutive equations in terms 
of the independent variables in the first row. 

2 Equations of thermoelasticity 

2.1 Constitutive relations of thermomechanics 

The thermomcchanical variables are the bulk stress and strain, cr and e, the temperature 
deviation, 6 from the ambient absolute temperature da (|^|/^a <C 1), and the entropy 
deviation per unit volume, s, from the ambient entropy s^. All quantities are defined relative 
to their ambient values, and would be zero in the absence of exterior motivating forces. The 
constitutive relations for strain and entropy in terms of the independent variables stress and 
temperature, {cr, 6*}, are [T71[T51[TU| 

e = Scr + ae', s^Cp{d/ea) + OL-rT. (3) 

Table ^ summarizes the alternative formulations of the same equations based on different 
choices of the independent variables: {e, 9}, {e, s] or {cr, s}. 

The material constants are as follows: S is the fourth order tensor of isothermal com- 
pliances, with inverse C, and corresponding adiabatic quantities are Ss and C^. The sym- 
metric second order tensor of thermal expansion coefhcients is a, and the related tensor 
/3 is defined as ,9 = Ca, while the quantities Cp and Cv are the heat capacities per unit 
volume at constant stress and strain respectively. The following relations can be verified 
from Table HI 

Cp/ea = {Cjea)+ofCcx, S = S, + (0„/Cp)a®a, C = C, -(0ja)/3®/3. (4) 

A word about notation: a ■ cr =tr(Q;cr) is a scalar, while a (g) a is a fourth order tensor. 
The constitutive relations in Tabled follow from the standard thermodynamic relations 

dU = (7 ■ de + 9ads, dF = cr ■ de — Sad9, (5) 



dG = e ■ dcr — Sad9, dH = — e • dcr + 9ads, 



(6) 
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where U, F, G and H are, respectively, internal energy, Helmholtz free energy, Gibbs free 
energy and enthalpy, all per unit volume. These are related by the standard connections 
U = F + TS = H + cre = G + TS + (t -e, where here, T and S are the absolute temperature 
and entropy, T = 0a, + O,S = Sa + s. The energy densities can be expressed, in the quadratic 
approximation that is used here, as 

2U = e-Ce+{Cjea)e^, 2F = e-Cse-{6a/C^)s^ (7) 



2G ^ -CT ■Ss(T ~ iea/Gp)s^, 2H ^-(T-S(T+{Gp/9a)9'^, (8) 

The constitutive relations in Table ^ follow from ((SJ and Q combined with the basic 
definitions of the thermal expansion coefficients and the heat capacity, 

de 



AQp,. = Cp,,„AT. (9) 



2.2 Thermoelastic relaxation governing equations 

We first present the exact governing equations and then make appropriate asymptotic ap- 
proximations. The motion is assumed to be caused by external forcing with no internal 
applied body forces or sources of heat. The heat flow in an elastic body is governed by the 
energy balance 

6l„.s + divq = 0, (10) 
where q is the heat flux. The equation for s in terms of 9 and cr in Table ^a-nd (|10|l imply 

Cp9 + diYCi = -9aOL- & (11) 

Irreversibility is introduced by requiring the heat flux to satisfy a generalized form of 
Fourier's relation [1^ 

q + r,.q + KVe = 0, (12) 

where K is the positive definite thermal conductivity tensor, and is the thermal relaxation 
time ^ni- The Cattaneo-Vernotte heat flux equation (|12|l includes the classical and more 
commonly used Kirchhoff law in the limit as Tr — > 0. The parameter is sometimes 
introduced to ensure finite speeds in the theory [2]]. Our objective is to solve the linear 
system of partial differential equations, (|ll|l and (|12|) . for 6* as a function of the forcing in 
the right member of Hll() . The heat flux can be eliminated to give a single equation for 9, 

Gp (9 + Tr9j - div KV9 = -(1 + Tr^)9a a ■ &. (13) 

A closed system of equations is obtained by applying the dynamic equilibrium condition 

diver - pii = 0, (14) 

where p the mass per unit volume and u is the elastic displacement vector, related to the 
strain via e = (Vu + (Vu)-^)/2. It is shown in Appendix A that closed-form solutions of 
the coupled system of equations H13|) and (I14|) are generally feasible only under restricted 
conditions. These require, essentially, that the material must be elastically isotropic, which 
is too limiting for our purposes. 
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3 Solution for arbitrary structures 
3.1 Asymptotic approximation 

The key quantities arc the positive definite compliance and stiffness differences AS = 
S — Sg and AC = Cg — C. Tliese determine tfie energy decrement between tfie final 
(isothermal) relaxed and initial unrelaxed (adiabatic) states, and they follow from Q as 
AS = {9a/Cp)a.oF and AC = Cg — C = {9a/Cy)^(3^ . Zener's approach is based on 
a separation between the mechanics and the thermodynamics. By assumption, the total 
difference between the relaxed and unrelaxed energies is small. Specifically AEq/Eq <^ 1, 
where the mechanical energy £q = ^ctq ■ eg is defined by eg and <Tg, which are related by 
the purely mechanical equation ctq — Ceg (ignoring temperature and entropy variations). 
Thus, £o = i(Tg • ScTg and the decrement may be defined as 

Afo = (Tg-AS<To = (^?a/Cp)(a-*To)'. (15) 

Alternatively, A£g w eg • ACeg = {0a/Cy){f3 ■ eg)^, where the approximation is due to the 
assumed purely mechanical relation between erg and eg. The main point is that the relative 
change in energy between the unrelaxed and relaxed states in either case is the same to 
leading order in e, where the nondimensional parameter governing TE damping is 

e = EdaayCp. (16) 

This definition of e is chosen to equal the relative change in elastic moduli Ac/cg for TE 
relaxation of a thin beam, equation It can also be expressed as 

e=i(l-2i.)(Cp-a)/Cp, (17) 

where v is the isothermal Poisson's ratio. Chadwick '15' employed a slightly different 
nondimensional parameter (denoted here as Ec to distinguish it from e) 

It is clear that the nondimensional parameters are closely related, and in particular, of the 
same order of magnitude. 



3.2 Solution by projections for anisotropic systems 

The coupled equations (|13|l and H14|l are solved using a regular perturbation procedure in 
the asymptotic parameter e ^ 1. We will achieve the solution using a projection method, 
similar to Zener's approach. A separation of variables reduces the problem to coupled 
ordinary differential equations in time. Anisotropy in the elastic material does not permit a 
modal expansion with a common set of scalar eigenfunctions, the basis for Zener's method, 
and the key to a generalization of his method to the limited but important case of isotropic 
solids see Appendix A. However, even in the case of the exact solution obtained by 
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Chadwick |15| , the interesting phenomena are obtained by the leading order approximation 
to the complex-valued frequencies. It therefore makes more sense ultimately to proceed by 
a regular asymptotic approximation at the stage of the coupled equations H13|l and H14|) . In 
this approach we view them as decoupled to leading order, whereby the elasticity problem 
is solved with no thermal effects. That is, we consider the elasticity as an uncoupled but 
vital forcing term in the 'thermal equation' H13() . 

Thus, to leading order equation H13|l gives an uncoupled equation for temperature with 
the stress entering on the right hand side as a forcing term. It is important to note that the 
forcing in H13|l is proportional to a • o". When the thermal expansion tensor is isotropic, the 
forcing depends upon the rate of hydrostatic stress, even when the material is elastically 
anisotropic. Thus, it is the hydrostatic stress, not strain, that governs the TE loss. 

Further progress is made using projections onto a set of eigenfunctions. In fact, this 
is similar to the method first proposed by Zener |Hj, which treated a simpler decoupled 
thermoelasticity problem. We first discuss the generalization of Zener 's method to H13|l as 
it allows us to determine the final answer in a form similar to the familiar and classical result 
of Zener for an elastic beam in flexure. We will later compare the general solution with 
direct solutions for particular configurations. Assume the temperature can be represented as 

oo 

0(x,i) = ^0„(t)(/)„(x), (19) 

where x = (z, y, z) = (xi, X2, ^3) and the eigenvalues t„ and eigenfunctions 0„ satisfy 

C-MivKV(^„+T-Vn = 0, 1, 2,..., (20) 

plus appropriate boundary conditions (for example, no flux). The amplitudes solve 

On + rJn + T-Hn = - (1 + ^-^) " • (21) 

where the brackets indicate the inner product (/, .g) = J dV f{x.) g{x.). 



3.3 A general result for energy loss 

Before considering applications of (|19|I - H21() to particular structures we first derive a general 
result for TE dissipation. A measure of local structural damping may be defined in terms 
of the local relative loss in energy per cycle. The rate of change of local mechanical energy 
per unit volume is £ = cr • e. We assume periodic oscillation for cr and e and determine the 
loss in the mechanical energy through the coupling to irreversible thermal process, 9, also 
periodic. Using the relation for e in terms of cr and 9 in Tabled gives 

£ = cr Sa + (7 ■ a9, (22) 

Taking the average over a cycle, and using // = where the overbar indicates a time 
average, implies that the local irreversible energy loss rate per unit volume is 



00 

£(x) = - 9{x,t) a ■ o-(x,t) = ^ <j)n{x) 9n{t)a ■ &{x,t). 

n=0 



(23) 
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The total power dissipated is thus 

oo 

dVSi^) = -Y,0n{t){K,OC-6r) 

^ — n 



E 

n=0 



(l-u;V,T„)2+a;2r2 



7f «.^)2 (24) 



where we have assumed periodic motion of period 27r/a; and used (|21|l to express the time 
harmonic temperature field in terms of the stress. 

The loss factor is then given by Q^^ sa tan (5 — A£/(27r£o), where A£ = lid 27r/(jj is the 
loss per cycle and Eq is the total energy of oscillation. We find 

- 7^ E n 2 ""X^ 2 2 " • ^)'- (25) 

This provides a general formula for the TE loss in terms of the inhomogeneous stress. The 
Q of a particular mode may be straightforwardly obtained by integrating H25() over the 
volume of the oscillator and dividing by the total energy. This equation may alternatively 
be expressed in terms of the inhomogeneous strain, however, wc find the stress formulation 
more convenient. 

Equation (|25|l is one of the main results of the paper, as it provides a means to compute 
TE dissipation given a solution in terms of the inhomogeneous stress. 

We remark on the summation in H25() . If the first term in the infinite sum is dominant, 
as is often the case jH|, the sum can be truncated after only one term (n = 0). This gives a 
result very similar to except that the simple Lorentzian in the latter is replaced by the 
generalized Lorentzian amplitude 

^(^^) = n 2 "^"^2 ^ 2 2 - (26) 

Of course, this reduces to the classical Lorentzian when — 0, which has a maximum 
as a function of uj when lot = 1 . It is worth describing the properties of this generalized 
Lorentzian, in particular how Tr influences the maximum. For every > 0, A has a single 
maximum at a unique value oi lo — lj* defined by 

u;*r = [(l-4r + 16r2)i/2_i + 2r]^^V(^\/6) where r = r^/r. (27) 

Furthermore, uj*t > 1 for a restricted range of Tr- Specifically, 1 < liJ*t < ^4/3 for 
< Tj. < 2r/3, with uj*t equal to unity at the two extremes (r^ = and — 2t/3) and 
uj*T — ■y/4/3 for Tr = r/4. Conversely, uj*t < 1 for > 2t/3. In particular, for relatively 
large 3> t, the maximum is at uj*{TTrY^'^ — 1. The value of A at the maximum, Amax, 
increases monotonically from Amax = 1/2 when there is zero thermal relaxation, = 0, to 

Amax ~ (Tr/r)^/^ for Tr > T. 
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4 Thermoelastically damped orthotropic thin plates 
4.1 Thin plate dynamics 

We consider plates that are orthotropic with axes of symmetry coincident with the coordi- 
nate axes. Assume, with no loss in generality, that the thermal expansion tensor is diagonal 
a =diag(ai, a2, as). By virtue of the thin plate configuration we can ignore the normal 
stress cTzz, and employing the Kirchhoff approximation for the deformation, we have 

a ■ a- = aiGxx + "20-^ = [(1 + V2i)Ei aiK^x + (1 + i^i2)i?2 o:2Kyy] , (28) 

i — i^l2l^21 

where k is the curvature tensor with components Kxx, i^yy and Kxy The curvature is related 
to the transverse deflection of the centre-plane, w{x, y), by 

-'-'iSH- 

The quantity Vij is the Poisson ratio for strain in the j-direction caused by stress in the 
z-direction, and the two Poisson's ratios satisfy i^2i^'i — i'i2-E'2- The instantaneous potential 
energy density per unit area of a plate in flexure is 

SpE = ^ I -; {EikI^ + E2kI + 2v2lElHxxl^yy) + 4/iK2 I , (30) 

Z i — I^12J^21 J 

where /i is the in-plane shear modulus and I = (z, z) — /12. 



4.2 Asymptotic solution by projections for thin plates 

Our purpose here is to obtain a general solution for the TE loss based on the assumption 
that the transverse diffusion of heat can be ignored. This assumption is examined in Section 
|S1 where it is shown that provided the assumptions of thin plate theory and classical heat 
transport are satisfied, the transverse heat flow gives rise only to small corrections to the 
TE loss. 

The governing equation for the non-equilibrium temperature field is obtained by ignoring 
the transverse heat flow terms in p3|l . 

C, [e + T,.e) - = - (^1 + r.|) 0, a . dr (31) 

where the dependence of the local temperature on position is governed solely by the x- 
dependence of the prescribed stress field a where now x — (x, y) — {xi,X2)- Also, is the 
through-thickness element of the thermal conductivity tensor, which in keeping with the 
general orthotropic formulation, is K = diag(i4ri, K2, K3). 

The analysis leading to H24I) for the power dissipated by TE effects can be repeated with 
the proviso that the eigenfunctions (/>„ of the heat equation are now functions only of z, 
the thickness direction, and inner products should in this case be interpreted accordingly 
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(/,<?) ~ J dzf{z) g{z). Hence the analogue to (|24|l refers not to the total power lost, but 
instead to the rate of energy loss per unit area at position x; 

The temperature modes of interest are antisymmetric about the mid-plane |Sj with 

</>„ = (2//i)i/2sin(2n+l)^, r„ = (2n + 1)-^ n = 0,1,2,.... (33) 

Temperature modes that are symmetric in z have zero coupling to flexural stress compo- 
nents, since they are antisymmetric, and also to any membrane stresses, which are constant 
across the thickness. It is evident from (|32|) and 128|l that the thermal loss in flexure depends 
upon the quantities [2| 

/„ = (0„, z) V(z, z) = 96/[(2n + (34) 
Combining (|2Hll, (EH), and (PSj), we obtain 

f(a;) = -£:D»ss(K(a;))— ^ V/n-r; ^ \2 , 22 (^5) 



n=0 



where E = {Ei+ E2) is the average Young's modulus, a = (ai +a2)/2, k is the curvature 
tensor, and the 'dissipation energy density' (per unit area) Eoiss is given by 



£Diss{l^{x)) 



Ea^ 



(1 + V2l)Eia2 Kxx + (1 + iyi2)E2a2 Kyy'' ^ 



1 - J/12J^21 



(36) 



Equation is a key result pertaining to TE dissipation in structures which can be mod- 
elled as thin plates. Unlike most previous results, the predicted dissipation is inhomoge- 
neous. 

The local TE dissipation depends on position via the dependence on the local curvature 
tensor k{x). This aspect may be explored by defining the quantity p{k) = EoissiK) / £q{k)^ 
which gives a measure of the local TE energy dissipation relative to the local deformation 
energy. Expressing the total energy as twice the average potential energy by virtue of the 
virial theorem we find, 

,^ _ [Ed?{l~ Vl2V2l)\ ^ ^ [{^ + V2l)Eiai + {'^ + Vl2)E2a2 K.yy\^ 

EikI^ + E2Kyy + 2l'2lEiKxxHyy + 4/i(l - I^12J^2l)K^y 

The parameter p simplifies for materials of cubic symmetry, such as silicon, with Ei = 
E2 — E = E , ai — a2 and 1/12 = 1^21 = i^, and hence 

p(«) ^ ( {^x. + ^yvY .33. 

\ 1 - Z// kIx + l^ly + '^Vfixxliyy + 4i? V(l " '^^)l^ly 
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For isotropic materials, E = 2/i(l 
p{k) = 



and p{k,) becomes 
(trK)2 



(tr/«)2 - 2(1 - iy)detK' 



In this case, p{k) depends upon the two principal invariants of the curvature: tiK = 



(39) 



and det k 

Ki 



"xy 



Let Ki and K2 be the two principal curvatures, satisfying 



K2 = trK and — det/«, then 

l + v 



p{k) 



1 



K2) 



(k1 + K2) - 2(1 - J/)kiK2 



Thus, 



0<p{n)<2/{l~v), 



(40) 

(41) 
K2, and p 



with p(/«) = at locations where the plate is locally saddle-shaped, ki = 
achieves its maximum value at points where it is locally spherical, ki = K2. 

The bounds H41(l also apply to materials of cubic anisotropy, and occur in the same 
circumstances as for the isotropic material, as can be verified from (|38|l . 



4.3 TE loss factors for flexural modes and waves 

The loss factor of a particular flexural mode may be predicted via H35|) once the displacement 
field w{x) of the mode is known. The curvature tensor is first evaluated via the standard 
relation H29|l . The total energy lost from the mode per cycle is then calculated by integrating 
H35|l over the volume of the plate and time averaging. In the most common situation, the 
displacement field will be evaluated in frequency space as a complex quantity, and the time 
average is obtained in the usual way as f{t)g{t) — ^Re [f{uj)g*{uj)], where / is the Fourier 
transform. Thus the TE dissipation for mode a is given by 



71=0 



J dA SoissiK-) 



(42) 



where £"00- = / dA {Eke + £pe) is the total energy of mode a. 

The modal energy £ocr may be computed as either twice the average kinetic energy of the 
system or twice the average potential energy of the system, whichever is more convenient. 
The kinetic energy is often preferable when using experimental data because the second 
derivatives appearing in the curvature tensor can be ill behaved in the presence of noise. 

The expression for the TE loss factor given above is closely related to the loss factor 
given by Zener for a simple beam in fiexure. Zener gave the result (see (|2l), 



Qzcncr ~ 



EOaa^ 



UJTq 



Cp 1 + 



(43) 



where we have included only the first term in the infinite sum. Our results may thus be 
interpreted as 

_ / dA£niss{i^) 



/cr 



Zcncr 



(44) 
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where the quantity was cahed the modal participation factor (MPF) in the isotropic case 
studied in [J. The MPF is evidently closely related to the quantity p{k) analyzed above. 

The MPF simplifies considerably for an isotropic medium. Using similar manipulations 
as employed above, we find 

P^ = ^;;^Jd^(^ (45) 

a result which differs from that given in ^ by a factor of (1 — i^)^^. The reason for the 
difference is that it was assumed in 0] that flexural energy would be dissipated at the rate 
predicted by Zener for a beam. However, a plate undergoes more compressive stress than a 
beam and thus dissipates more energy accordingly. 

It is not difficult to show, based on 145() . that the MPF for travelling flexural waves in 
an isotropic plate is 

Pa = Y~ — for a flexural wave. (46) 

This shows the difference between the TE dissipation of a cantilever beam (p^ — 1) versus 
the corresponding vibration of a plate. The identity 14t)|) may be obtained by explicit 
substitution of a flexural waves solution, or more simply, as follows. Integrating by parts 
and using the governing plate equation EI{1 — i'^)^^W^w ~ phLo^w = 0, gives 

j dAjt^= J dAlJ^^^^^j^ J dAl^^^^^^j^So- (47) 

The last equality employs the expression for the potential energy, (|30() . Also, it is assumed 
that the plate boundary conditions may be ignored in the above integration, which is true 
for a travelling wave in a plate of 'infinite' extent. 



5 Effective thin plate equations and flexural waves 

The general theory is now applied to thin plates in flexure with the goal of deriving general 
governing equations similar to the classical Kirchhoff thin plate equations. Our objective 
is to provide an alternative means of calculating p^ for flexural waves, and also to examine 
the range of validity of our approximation in ignoring lateral thermal diffusion. 
The equation for 6, H13|l . becomes for time harmonic motion (e~"^* assumed) 

where K = diag {Ki, K2, K3) and 

k'^ = iuj{l - iLUTr)Cp/K3. (49) 

If the stress term on the right hand side of H48I) is weakly dependent on x and y, then we may 
argue that 9 inherits the same weak dependence. Provided this dependence on transverse 
position is uniform, the case of simple plane wave, the final term may be combined with the 
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k^9 term, and hence interpreted as modifying the thermal diffusion rate. Here, we ignore 
this and simphfy H48|l to 



(50) 



The importance of the simphfication (|50|l at this stage is that it allows us to derive a set of 
effective plate equations in which the TE damping appears directly. This approach follows 
on that of Alblas J.2 and of Lifshitz and Roukes who derived the effective equation 
governing the motion of a thermoelastically damped beam. 

Assuming zero flux conditions a,t z = ±/i/2, and noting that the right member of (|50|l 
is proportional to z in flexure, cr = (z, (t)I~^z, we find that the solution is 



/a 



(z, a ■ cr) 



sin kz \ 
kcos{kh/2) J 



(51) 



The thermally perturbed stress is therefore, using the expression for entropy in terms of 
temperature and strain from Tabled 



Ce 



/a 



(Z, OL ■ (To) 



sin kz \ 



k coa{kh/2) ^ 



(52) 



The stress in the absence of thermal effects, ctq, is proportional to z. In order to apply 
l)52|l to the thin plate the standard plane stress conditions must be enforced. We consider 
an orthotropic plate with a symmetry plane coincident with the neutral plane (z = 0), for 
which the standard stress/strain relations for plane-stress are 



XX 



El 



1 — 1/121^21 1 — 1'12'/21 
^12-^2 E2 



1 — 1/12!/21 1 — I^12f21 

2^ 









^XX 








^yy 











(53) 



where a^J , ei'i^ etc. are the stresses and strains in the absence of the thermoelastic damp- 
ing. Based on (|52|l this implies that the in-plane TE stresses are 



.(0) 



^yy — "yy 



>xy — f^xy ■ 



(1 - iyi2:^2i)CpI ^ ' ^ +(^2<Jyy)^z ^^^^ 



sin kz \ 
Jkhj2))' 



{a2 + aiVi2)eaE2 1 (0) , (o)\ / sinfcz \ 

[l — i^i2i^2i)CpI \ kcos[kh/2) J 



This provides the variation of stress through the thickness due to the temperature varia- 
tion. The moments are found by taking the first moment of the stresses through the plate 
thickness, leading to 



Mx 



(q1 + a2'^2l)0aEi 
(1 - Vl2V2l)Cp 



f(kh) (aiMW+a2M(°)), 
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^ M(?) + + "^"^^y («iMW + a.Af W), (55) 



r(o) _ ^(0) 



yy 

xy 1 



(1 - i^l2!^2l)C'p 



where M^ia;-^ = (z, CTxa; ), etc., and the function / is 



/(C) = 1 + 7J 



— — tan — 
2 2 



(56) 



The standard moment-curvature relations follow from ()53fl and the Kirchhoff assumption 
as 









= / 







l-iyi2;^21 1-I'12l'21 

2^l 

The governing equation for the thin plate is 



-El 



1-1^121^21 l-I'12l'21 



1^12-^2 



i?2 













r 1 




Ky 







(57) 



+ phu) 1/7 = 0. 



dx^ dxdy dy"^ 

Thus, we obtain the effective thin plate equation for w 

--. [EiWxxxx + E2Wyyyy + 2[2^(l - i^l2J^2l) + '^2lEi]wa;. 

i — J^12i^21 L 

fikh)9al 



(58) 



xyy 



(1 - i^l2i^2l)2Cp 



(ai + a2i^2i) E^Wxxxx + 2(ai + a2:^2i)(a2 + aiiyi2)EiE2'w^ 



+ + aii^l2) E2Wyyyy > + pkUJ W = 



For a plate made of material with cubic symmetry, this reduces to 



xxyy 



(59) 



f{kh) 



^1^ + 4(1 - v) 



^ 'e 2 



"xxyy —=-—phuj W = 0. 



(60) 

The second term vanishes for isotropic materials, in which case we may combine the damping 
with the plate stiffness to get an equation in the standard form, 

DV^w ~ phuj"^ w ^ 0, (61) 

where the TE loss is now contained in the complex-valued flexural stiffness 



D = 



EI 



1 



f{kh) 



(62) 
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In the more general orthotropic situation, the TE damping is inhomogeneous, as in H59() . 
and cannot be interpreted in terms of a single frequency-dependent complex modulus of 
elasticity. 

Our equation (|51|l . a generalization of Zener's results for a beam p , gives a homogeneous 
damping as a result of assuming a uniform plane wave vibration. This provides some 
guidance as to the loss factors of vibrating thin plate structures but finite structures will 
support resonant modes that consist of a variety of wavenumbers which can interfere with 
one another, and thus the actual loss factor for a particular mode may differ significantly 
from this value. A particularly salient example is the case of a twisting mode for which 
the source field for temperature fluctuations, trir (for isotropic cr), is very small, and the 
resulting values of may be orders of magnitude smaller than the result given in H61|l . 

Corresponding results for flexural waves in circular rods are in Appendix C. In particular 
we note that the function / takes on a different form from that for a plate. 

5.1 Dispersion relation including in-plane variation 

It is now relatively straightforward to revise the analysis of Section ^ to include in-plane 
variation in both the stress and the temperature. We begin by assuming that all field 
variables possess in-plane dependence e*""^ with wavenumber C. Then H48() becomes 



where 7 = (fc^ — (^^ Ki/ K^y^^ and k is defined in (|49|l . Solving as before, we find that the 
temperature is 



and similar generalizations can be obtained for the stresses and moments, H52ll - H57(l . 

We cannot derive a governing equation, similar to for example, since now the 

in-plane dependence of the stress and hence w has been assumed a priori. However, the 
in-plane wavenumber can be obtained in a self-consistent manner from the latter equation 
by assuming w — wqc*^^, where wq is constant. This yields an equation for ^, 



(63) 




(64) 



1 + z/\ E0aa^ 
l^J Cp ^ 



,2 . -1 



(65) 



the solution of which we discuss next. 



5.2 Effects of transverse TE dissipation 

The wavenumber of a flexural wave in an undamped thin plate is fc/, where 

kj =uj^{l-v^)ph/{EI), 



(66) 



Treating C as an asymptotic series in the small parameter e, it is clear that the leading order 
solution to the general dispersion relation of 165() is C = ^/+0(£)- The next term is given 
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by 



e 

1 - - 



7/ 



4 V 1 - J/ 



(67) 



where 7/ = (fc^ — k'jy^'^. Direct substitution gives 

7/ = /co (1 + «a^mfp//i - it^T^)^^^ , (68) 
where /cq = {—iioCp/ K^Y^'^ and ?mfp; the mean free path for phonons at temperature T, is 

3i^3(T) 



'mfp(r) — 



(69) 



The quantity c is the average elastic wave speed, and the order one constant a is given by 

a = [2(1 - iy)p/{3^)f^ (Ki/K^) c. (70) 

Values of the mean free path at room temperature are typically on the order of tens of 
nanometres using the value c = {cl + 2ct)/3 where ct — \/l^/p is the transverse wave 
speed and cl — [2(1 — — 2h')]^^'^ ct the longitudinal speed. For isotropic conductivity 
the parameter a is a function of Poisson's ratio, 



Vl — 2i/ 



/V27 



(71) 



which is of order unity. 

Equation (|68() provides an avenue to compute corrections to classical TE dissipation 
arising both from transverse diffusion and thermal relaxation. The asymptotic solution 
H67|l indicates a propagating flexural wave has attenuation equal to Im^. Using the relation 
between attenuation and Q, and noting that the undamped flexural wave has group velocity 
2uj/kf, we find 1/Q — Im4C/fc/, or 



l + v 



Im 



1 — iujTr 



1 — ILOTr + ialmlp/h 



(72) 



In order to proceed further in understanding this we replace the function / by its series 
expansion. Equations IjB.ip . (|B.2|) and IjB.Sp imply 



00 



1 



where t„, n = 0, 1, 2, . . . are defined in (|33|l 9. Thus, 



(73) 



l + v 
1 - V 



Tn (1 +UJTn{aln^ip/h)) 



n=0 



(74) 
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The classical result for the TE loss in a vibrating beam is |S] 



Q-'=^T^- (75) 

1 + UJ'^Tq 



This ignores the effect of t^, and for the purpose of comparison with H74|) we consider the 
case Tr = in the latter. We truncate H74() at the first term (which is reasonable considering 
/o//i = 81), and for further simplicity, the factor /o = 96/7t'^ — 0.9855 is replaced by unity, 
giving 

Comparing (|75|) and H7t)|l the first difference is the factor {1 + i/) / {1 — i/) which can be 
attributed to the fact that it is a plate rather than a beam as discussed in the previous 
subsection, see (|44f) and (|46|l . The major additional distinction between the present theory 
and the classical result is the presence of the terms involving Zmfp resulting from transverse 
diffusion, which are absent in Zener's theory and previous analyses. These corrections are 
always small, since the domain of validity of the thermodynamic analysis employed here is 
restricted to distances far greater than the mean free path, the domain h 3> Zmfp- 

Alternatively, we note that 176() contains two characteristic times: the 'Zener' relaxation 
time To ~ h'^Cp/n'^K^, and a new characteristic time defined by r* ~ ''o(aZ,nfp//i). It 
is easily checked that ujt* — [kfh/ir)'^, which must be a small quantity in order for the 
Kirchhoff assumption to remain valid, that is, a necessary prerequisite for the validity of 
the plate theory is that the flexural wavelength is much longer than the thickness. Also, 



r* h7r-'^y/l2{l-i'^)p/E, (77) 

which suggests interpreting r* as the travel time of an elastic wave across the thickness. 
But this time must be much less than 1/w, otherwise thickness resonances can occur, again 
violating the conditions of the plate theory. 

Thus, the restrictions on the use of both the thermal conduction and the thin plate 
theories requires that u!TQ{aln-ifp/h) is small. Nevertheless, small corrections arising from 
transverse diffusion can be estimated in this fashion. 

The effect of non-zero can be considered by using an estimate for this relaxation time. 
Rudgers |2J provides the estimates = Tmfp/3, where Tmfp = Imip/c (Rudgers actually 
calculates two approximations for but they are of the same order of magnitude). The term 
ujTr must remain small in the context of thin plate theory, otherwise the same assumptions 
as before are violated, for example, the wavelength is far less than the thickness. However, if 
the phonon mean free path is comparable with the plate thickness, then both the transverse 
diffusion and the Cattaneo-Vernotte terms can become important. 



6 Conclusion 



Equation H25(l is one of the basic results of the paper, as it provides a means to compute TE 
dissipation given a solution in terms of the inhomogeneous stress. This computation is in 
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general a complicated undertaking because it requires the determination of the cigcnfunc- 
tions and eigenvalues of the heat equation in the particular geometry of interest. Never- 
theless, a clear recipe for an approximate calculation of the TE dissipation in an arbitrary, 
anisotropic, elastic medium is given. 

The theoretical scope is then narrowed to the case of thin plate structures for which the 
thermal heat flow simplifies dramatically. In this case, we are able to obtain an explicit 
formula for the TE loss in an arbitrary, anisotropic elastic system, (|36ll . This equation, which 
shows that TE loss is in general inhomogeneous, is the other principal result obtained. It 
should be emphasized that the A£ derived in 1)36(1 has the property that it is not a global 
measure of damping but is local, and can be used to define TE loss at a point in a structure, 
and then integrated to determine the Q of an arbitrary vibrational mode. This distinction 
is particularly important for geometries of interest in MEMS/NEMS applications, for the 
TE loss rate of low order modes in typical geometries may vary significantly with position. 
The local result for TE loss predicts that the attenuation of a flexural wave in a large plate 
is (1 + v) / {1 — v) times the attenuation of a wave in a beam at the same frequency. 

Finally, we have investigated the effects of transverse diffusion. Transverse diffusion 
effects are most easily analyzed for plane wave behaviour. For this case we have derived an 
effective plate equation including the effects of thermoelasticity to leading order in the TE 
coupling. The TE dissipation can be examined fairly simply in the plane wave case because 
TE losses are homogeneous, and are therefore contained in a single loss factor. Corrections 
to TE loss resulting from transverse diffusion arc found to be generally small within the 
domain of validity of our theoretical models. 
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A General solutions of the equations of thermo elastic- 
ity 

The governing equations and 114(1 are more commonly expressed as coupled equations 
of temperature and displacement. This can be achieved by first eliminating cr explicitly: 

divCe - pii - /3V6I = 0, (A.l) 

e + Tre-c-^divKye + {i + Tr^)^f3-e = o. (a.2) 

These may be expressed more succinctly as follows for time harmonic motion (e~''^* is 
assumed) , 

C{V , Lu)lJ + Lu^V = 0, (A.3) 

where 



p-1 Q(V) iLuipCJOa)-'/^ b(V) 

ILuipCJea)-'/^ b^(V) ^(V, OU) 



(A.5) 



Qiji^) = Gikji VkVi, b(v) = /3v and k{v,uj) = [9a{Tr + {—icLj)~^)]~^ v • Kv. Equation ((A.3|l 
represents a generalized eigenvalue problem for the complex- valued modal frequencies. 

Some simplification is possible for isotropic bodies, for which a = al and (3 = SktoI 
where kt = E/{1 — 2v) is the isothermal bulk modulus. Assuming the ansatz |15| 

U = j Tf 1 , ^here V^i) + A^V^ = 0, (A.6) 



the eigenvalue problem l|A.3l) becomes 



2 2 a2 I iippX 



2 



^ -TTT^^— wTt M = 0, (A.; 



where cl is the longitudinal wave speed. Eliminating the constant A gives an equation for 
oj in terms of the wavenumber A: 

(-^ - ciA^) f ^.f^' . - zc) + .cA^^ = 0. (A.9) 

Thus, the problem of determining the eigenfrequencies is reduced to finding solutions to the 
Helmholtz equation l)A.6() 9 in the domain of interest. Chadwick ^Hl discusses this further 
and provides formulas for the roots of l|A.9l) for = 0, in which case it reduces to a cubic 
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The separation of variables approach does not work for generaUy anisotropic bodies. 
However, it is worth noting that solutions of the form (|A.6|I are valid as long as the thermal 
properties are isotropic and the elasticity tensor C satisfies Cij^ifijUkrii = Con^rii for all n, 
and some constant Cq . The most general form of anisotropy with this property has stiffness 
(using Voigt notation) 
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B Exact results for the function / 

Some exact results are presented for the function / of (|56|l occurring in the TE theory for 
beams and plates in flexure. We begin with the representation 



/(o = i+E- 



,=0 2^^^-l 



(B.l) 



(B.2) 



The infinite series in IjB.ljl is a consequence of the fact that the left member is a meromorphic 
function of ^ with simple poles at ^ = ±(2j + 1)tt, and residues that are readily determined. 
Note that 



where fj, j = 0, 1, 2, . . ., are defined in and 



It follows from ljB.l|) that for real- valued ^, 



(B.3) 



We also note the identity 



r,,. , -, 6 r (l~t) ^ sinh^-isin^ ^ ^ 

+ = '-eV^ -^^ coshg + cosg ) ^ ' 



(B.4) 



(B.5) 
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which is useful in the case of — (r — 0), as it provides a closed form expression for the 
TE damping in that case via |13| 

,0.-l{i-l(|^!li±|Li)}. (B.) 

The exact result for a beam in flexure was derived by Lifshitz and Roukes |13| . and 
corresponds to the use of the function / of H56I) . The equivalence of their derivation and 
Zener's prediction is confirmed by (|B.4|I (both Zener and Lifshitz and Roukes considered 
the case r = 0). Thus, while the closed form expression of Lifshitz and Roukes is novel, it is 
simply a more concise expression of the infinite series of Zener. Both are based on the same 
implicit or explicit first order approximation arising from the small difference between the 
adiabatic and isothermal systems. 

In fact, the analysis of Lifshitz and Roukes TS' is a special case of the more general 
problem solved by Alblas |jl2|, in which the lateral (third) dimension of the beam is taken 
into account and the boundary value problem for the temperature is more general. How- 
ever, Alblas considered the situation in which the boundary condition for the temperature 
variation is zero, and as a result he obtained a different functional form of the TE damping. 



C Thermoelasticity theory for a circular rod 



In the case of a circular rod of radius a we find that the solution to H5U|I is 



(C.l) 



where Ic = {z, z) — z = rcos?/j and fc is defined in H49|l . The function associated 

with the moment of the temperature is f{ka) where now 



/(C) = 1 + 



c 



MO 



(C.2) 



Zener |22| analysed the circular rod using the projection method involving a series rather 
a closed form expression. By comparing the results here with those of Zener, we conclude 
that Eg 



fn 



iliil - 1)' 



where J( (<?„) =0, n = 0, 1, 2, . 
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